{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import time\n",
    "import datetime\n",
    "start = time.time()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import scanpy as sc\n",
    "import matplotlib.pyplot as plt\n",
    "import glob\n",
    "import seaborn as sns\n",
    "import bbknn\n",
    "from pandas.core.index import RangeIndex\n",
    "from matplotlib import rcParams\n",
    "from matplotlib.pyplot import rc_context"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.settings.verbosity = 3  \n",
    "sc.logging.print_header() \n",
    "sc.settings.set_figure_params(dpi=80)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V2IA = sc.read_loom('./V2h_Amp.loom')\n",
    "V2IU = sc.read_loom('./V2h_Isth.loom')\n",
    "P2IA = sc.read_loom('./P2h_Amp.loom')\n",
    "P2IU = sc.read_loom('./P2h_Isth.loom')\n",
    "V24IA = sc.read_loom('./V24Am.loom')\n",
    "V24IU = sc.read_loom('./V24Is.loom')\n",
    "P24IA = sc.read_loom('./P24h_Amp.loom')\n",
    "P24IU = sc.read_loom('./P24h_Isth.loom')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V2IA.var_names_make_unique()\n",
    "V2IU.var_names_make_unique()\n",
    "P2IA.var_names_make_unique()\n",
    "P2IU.var_names_make_unique()\n",
    "V24IA.var_names_make_unique()\n",
    "V24IU.var_names_make_unique()\n",
    "P24IA.var_names_make_unique()\n",
    "P24IU.var_names_make_unique()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V2IA.write('./V2hInfAmp.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V2IU.write('./V2hIsthUTJ.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P2IA.write('./P2hInfAmp.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P2IU.write('./P2hIsthUTJ.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V24IA.write('./V24hInfAmp.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V24IU.write('./V24hIsthUTJ.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P24IA.write('./P24hInfAmp.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P24IU.write('./P24hIsthUTJ.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "V2IAIU = V2IA.concatenate(V2IU, batch_key='Region', batch_categories=['InfAmp','IsthUTJ'])\n",
    "P2IAIU = P2IA.concatenate(P2IU, batch_key='Region', batch_categories=['InfAmp','IsthUTJ'])\n",
    "V24IAIU = V24IA.concatenate(V24IU, batch_key='Region', batch_categories=['InfAmp','IsthUTJ'])\n",
    "P24IAIU = P24IA.concatenate(P24IU, batch_key='Region', batch_categories=['InfAmp','IsthUTJ'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2hIAIU = V2IAIU.concatenate(P2IAIU, batch_key='Treatment', batch_categories=['V','P'])\n",
    "VP24hIAIU = V24IAIU.concatenate(P24IAIU, batch_key='Treatment', batch_categories=['V','P'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU = VP2hIAIU.concatenate(VP24hIAIU, batch_key='Time', batch_categories=['2h','24h'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.write('./VP2h24hIAIU.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.settings.verbosity = \n",
    "sc.logging.print_header() \n",
    "sc.settings.set_figure_params(dpi=80, facecolor='white')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pwd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU = sc.read( './VP2h24hIAIU.h5ad') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.var_names_make_unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.pl.highest_expr_genes(VP2h24hIAIU, n_top=20)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.filter_cells(VP2h24hIAIU, min_genes=200) \n",
    "sc.pp.filter_genes(VP2h24hIAIU, min_cells=3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "mito_genes = VP2h24hIAIU.var_names.str.startswith('mt') \n",
    "VP2h24hIAIU.obs['pct_counts_mt'] = np.sum(\n",
    "    VP2h24hIAIU[:, mito_genes].X, axis=1).A1 / np.sum(VP2h24hIAIU.X, axis=1).A1 \n",
    "VP2h24hIAIU.obs['total_counts'] = VP2h24hIAIU.X.sum(axis=1).A1 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.pl.violin(VP2h24hIAIU, ['n_genes'], jitter=0.4)\n",
    "sc.pl.violin(VP2h24hIAIU, ['total_counts'], jitter=0.4)\n",
    "sc.pl.violin(VP2h24hIAIU, ['pct_counts_mt'], jitter=0.4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.scatter(VP2h24hIAIU, x='total_counts', y='pct_counts_mt')\n",
    "# sc.pl.scatter(VP2h24hIAIU, x='total_counts', y='n_genes_by_counts')\n",
    "sc.pl.scatter(VP2h24hIAIU, x='total_counts', y='n_genes')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU = VP2h24hIAIU[VP2h24hIAIU.obs.n_genes < 4000, :]\n",
    "VP2h24hIAIU = VP2h24hIAIU[VP2h24hIAIU.obs.pct_counts_mt < 0.1, :]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.normalize_total(VP2h24hIAIU, target_sum=None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.log1p(VP2h24hIAIU)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.regress_out(VP2h24hIAIU, ['total_counts', 'pct_counts_mt'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.highly_variable_genes(VP2h24hIAIU, min_mean=0.0125, max_mean=3, min_disp=0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.highly_variable_genes(VP2h24hIAIU)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.raw = VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU = VP2h24hIAIU[:, VP2h24hIAIU.var.highly_variable]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.write('./CheckPoint.h5ad')\n",
    "CheckPoint = sc.read('./CheckPoint.h5ad')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.scale(VP2h24hIAIU, max_value=10) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.tl.pca(VP2h24hIAIU, svd_solver='arpack')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.pca(VP2h24hIAIU, color='Treatment')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.pl.pca(VP2h24hIAIU, color='Time')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.pca_variance_ratio(VP2h24hIAIU, log=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.pca(VP2h24hIAIU, wspace=0.35,\n",
    "          color=['Treatment', 'Time', 'Region', 'Epcam', 'Ovgp1', 'Foxj1', 'Myh11', 'Cd52', 'Pecam1', \n",
    "                  \n",
    "                  'Pdgfra', 'Twist2', 'Igf1', 'Dcn', 'Serpina1e', 'Myb', 'Pax8', 'Pax2', 'Cxcl17', 'Id2', 'Pdxk', 'Ccdc153'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.pca_variance_ratio(VP2h24hIAIU, n_pcs=40, log=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "\n",
    "sc.pp.neighbors(VP2h24hIAIU, n_neighbors=10, n_pcs=40)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.tl.leiden(VP2h24hIAIU)\n",
    "sc.tl.paga(VP2h24hIAIU)\n",
    "sc.pl.paga(VP2h24hIAIU, plot=False)  # remove `plot=False` if you want to see the coarse-grained graph\n",
    "sc.tl.umap(VP2h24hIAIU, init_pos='paga')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.umap(VP2h24hIAIU)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.tl.leiden(VP2h24hIAIU, resolution=0.2) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU,color=['leiden','Treatment', 'Time', 'Region'],legend_loc='on data')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU,\n",
    "           color=['leiden', 'Treatment', 'Time', 'Region', 'Epcam', 'Ovgp1', 'Pax8', 'Pax2', 'Foxj1', 'Ccdc153', \n",
    "                  'Myh11', 'Cd52', 'Ptprc', 'Pecam1', \n",
    "                  \n",
    "                  'Pdgfra', 'Twist2', 'Igf1', 'Dcn', 'Serpina1e', 'Myb', 'Cxcl17', 'Id2', 'Pdxk'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.obs.groupby([\"leiden\",\"Treatment\",\"Time\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.obs.groupby([\"leiden\",\"Treatment\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.obs.groupby([\"leiden\",\"Time\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU.obs.groupby([\"leiden\",\"Region\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "   \n",
    "VP2h24hIAIU_batchcorrected = bbknn.bbknn(VP2h24hIAIU, batch_key='Time', copy=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "sc.pl.pca(VP2h24hIAIU_batchcorrected, wspace=0.35,\n",
    "          color=['leiden', 'Treatment', 'Time', 'Region', 'Epcam', 'Ovgp1', 'Foxj1', 'Myh11', 'Cd52', 'Pecam1', \n",
    "                 \n",
    "                  'Pdgfra', 'Twist2', 'Igf1', 'Dcn', 'Serpina1e', 'Myb', 'Pax8', 'Pax2', 'Cxcl17', 'Id2', 'Pdxk', 'Ccdc153'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "sc.tl.umap(VP2h24hIAIU_batchcorrected)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, wspace=0.35, \n",
    "           color=['leiden', 'Region', 'Time', 'Treatment', 'Epcam', 'Ovgp1', 'Foxj1', 'Myh11', 'Cd52', 'Pecam1', \n",
    "                  \n",
    "                  'Pdgfra', 'Twist2', 'Igf1', 'Dcn', 'Serpina1e', 'Myb', 'Pax8', 'Pax2', 'Cxcl17', 'Id2', 'Pdxk', \n",
    "                  'Ccdc153', 'Ptprc'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "marker_genes = ['Epcam', 'Krt8','Krt19','Vim','Pdgfra','Twist2', 'Sprr2f',\n",
    "                'Acta2', 'Myh11','Foxj1', 'Ccdc153', 'Ovgp1','Wt1','Serpina1e',\n",
    "                'S100g','Pecam1', 'Cd52', 'Pgr', 'Esr1', 'Crabp2', 'Cd44', 'Pdxk'\n",
    "                 ]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color=marker_genes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.obs.groupby([\"leiden\",\"Treatment\",\"Time\",\"Region\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.obs.groupby([\"leiden\",\"Treatment\"]).apply(len)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_cluster_names = [\n",
    "     '0 Secretory_InfAmp',\n",
    "     '1 Secretory_IsthUTJ',\n",
    "     '2 Ciliated_InfAmp',\n",
    "     '3 Fibroblast_Pdgfra+',\n",
    "     '4 Fibroblast_Pdgfra-',\n",
    "     '5 Secretory+Ciliated_Veh',\n",
    "     '6 Secretory_IsthUTJ_subset1',\n",
    "     '7 Immune1',\n",
    "     '8 Cliated_IsthUTJ1',\n",
    "     '9 Cliated_IsthUTJ2',\n",
    "     '10 Endothelial1',\n",
    "     '11 Secretory_IsthUTJ_subset2',\n",
    "     '12 Muscle',\n",
    "     '13 Immune2',\n",
    "    '14 Endothelial2'\n",
    "]\n",
    "VP2h24hIAIU_batchcorrected.rename_categories('leiden', new_cluster_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(8, 2))\n",
    "for i in range(28):\n",
    "    plt.scatter(i, 1, c=sc.pl.palettes.zeileis_28[i], s=200)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "zeileis_colors = np.array(sc.pl.palettes.zeileis_28)\n",
    "new_colors = np.array(VP2h24hIAIU_batchcorrected.uns['leiden_colors'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_colors[[0]] = zeileis_colors[[0]]  \n",
    "new_colors[[1]] = zeileis_colors[[7]]  \n",
    "new_colors[[2]] = zeileis_colors[[5]]\n",
    "new_colors[[3]] = zeileis_colors[[16]]  \n",
    "new_colors[[4]] = zeileis_colors[[17]] \n",
    "new_colors[[5]] = zeileis_colors[[6]] \n",
    "new_colors[[6]] = zeileis_colors[[18]] \n",
    "new_colors[[7]] = zeileis_colors[[12]] \n",
    "new_colors[[8]] = zeileis_colors[[4]] \n",
    "new_colors[[9]] = zeileis_colors[[3]] \n",
    "new_colors[[10]] = zeileis_colors[[26]] \n",
    "new_colors[[11]] = zeileis_colors[[19]] \n",
    "new_colors[[12]] = zeileis_colors[[24]] \n",
    "new_colors[[13]] = zeileis_colors[[13]] \n",
    "new_colors[[14]] = zeileis_colors[[27]] "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.uns['leiden_colors'] = new_colors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color='leiden', legend_loc='right margin', legend_fontsize='8.0',\n",
    "           frameon=False, save='VP2h24hIAIU_AllPops_Leiden_060724.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color='Time', legend_loc='right margin', legend_fontsize='8.0',\n",
    "           frameon=False, use_raw=False, save='VP2h24hIAIU_AllPops_Time_060724.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_cluster_names = [\n",
    "    '2h',\n",
    "    '24h'\n",
    "]\n",
    "VP2h24hIAIU_batchcorrected.rename_categories('Time', new_cluster_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(8, 2))\n",
    "for i in range(28):\n",
    "    plt.scatter(i, 1, c=sc.pl.palettes.zeileis_28[i], s=200)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "zeileis_colors = np.array(sc.pl.palettes.zeileis_28)\n",
    "new_colors = np.array(VP2h24hIAIU_batchcorrected.uns['Time_colors'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_colors[[0]] = zeileis_colors[[6]]  #0 2h\n",
    "new_colors[[1]] = zeileis_colors[[17]]  #1 24h"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.uns['Time_colors'] = new_colors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color='Time', legend_loc='right margin', legend_fontsize='8.0',\n",
    "           frameon=False, save='VP2h24hIAIU_AllPops_Time_060724.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_cluster_names = [\n",
    "    'InfAmp',\n",
    "    'IsthUTJ'\n",
    "]\n",
    "VP2h24hIAIU_batchcorrected.rename_categories('Region', new_cluster_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(8, 2))\n",
    "for i in range(28):\n",
    "    plt.scatter(i, 1, c=sc.pl.palettes.zeileis_28[i], s=200)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "zeileis_colors = np.array(sc.pl.palettes.zeileis_28)\n",
    "new_colors = np.array(VP2h24hIAIU_batchcorrected.uns['Region_colors'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_colors[[0]] = zeileis_colors[[9]]  #0 InfAmp\n",
    "new_colors[[1]] = zeileis_colors[1]  #1 IsthUTJ"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.uns['Region_colors'] = new_colors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color='Region', legend_loc='right margin', legend_fontsize='8.0',\n",
    "           frameon=False, save='VP2h24hIAIU_AllPops_Region_060724.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_cluster_names = [\n",
    "    'Veh',\n",
    "    'P4'\n",
    "]\n",
    "VP2h24hIAIU_batchcorrected.rename_categories('Treatment', new_cluster_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(8, 2))\n",
    "for i in range(28):\n",
    "    plt.scatter(i, 1, c=sc.pl.palettes.zeileis_28[i], s=200)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "zeileis_colors = np.array(sc.pl.palettes.zeileis_28)\n",
    "new_colors = np.array(VP2h24hIAIU_batchcorrected.uns['Treatment_colors'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "new_colors[[0]] = zeileis_colors[[2]]  #0 Veh\n",
    "new_colors[[1]] = zeileis_colors[5]  #1 P4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected.uns['Treatment_colors'] = new_colors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected, color='Treatment', legend_loc='right margin', legend_fontsize='8.0',\n",
    "           frameon=False, save='VP2h24hIAIU_AllPops_Treatment_060724.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = sc.pl.dotplot(VP2h24hIAIU_batchcorrected, \n",
    "                   marker_genes, groupby='leiden', dendrogram=True, color_map='Blues',\n",
    "                  save='VP2h24hIAIU_AllPops_MarkerDotPlot_060724')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = sc.pl.dotplot(VP2h24hIAIU_batchcorrected, \n",
    "                   marker_genes, groupby='Region', dendrogram=True, color_map='Blues',\n",
    "                  save='VP2h24hIAIU_AllPops_MarkerDotPlot_Region_060724')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = sc.pl.dotplot(VP2h24hIAIU_batchcorrected, \n",
    "                   marker_genes, groupby='Treatment', dendrogram=True, color_map='Blues',\n",
    "                  save='VP2h24hIAIU_AllPops_MarkerDotPlot_Treatment_060724')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.rank_genes_groups(VP2h24hIAIU_batchcorrected, 'Treatment', method='wilcoxon',n_genes=1000)\n",
    "sc.pl.rank_genes_groups(VP2h24hIAIU_batchcorrected, n_genes=1000, sharey=False)\n",
    "df = pd.DataFrame(VP2h24hIAIU_batchcorrected.uns['rank_genes_groups']['names']).head(1000)\n",
    "df.to_csv('./VP2h24hIAIU_AllPops_Treatment_top1000.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "VP2h24hIAIU_batchcorrected_Epithelial = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['0 Secretory_InfAmp',\n",
    "                                                                                                                  '1 Secretory_IsthUTJ',\n",
    "                                                                                                                  '2 Ciliated_InfAmp',\n",
    "                                                                                                                  '5 Secretory+Ciliated_Veh',\n",
    "                                                                                                                  '6 Secretory_IsthUTJ_Subset1',\n",
    "                                                                                                                  '8 Ciliated_IsthUTJ1'\n",
    "                                                                                                                 '9 Ciliated_IsthUTJ2',\n",
    "                                                                                                                 '11 Secretory_IsthUTJ_subset2'])]\n",
    "\n",
    "VP2h24hIAIU_batchcorrected_Secretory = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['0 Secretory_InfAmp',\n",
    "                                                                                                                  '1 Secretory_IsthUTJ',\n",
    "                                                                                                                  '5 Secretory+Ciliated_Veh',\n",
    "                                                                                                                  '6 Secretory_IsthUTJ_Subset1',\n",
    "                                                                                                                 '11 Secretory_IsthUTJ_subset2'])]\n",
    "\n",
    "VP2h24hIAIU_batchcorrected_Ciliated = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['2 Ciliated_InfAmp',\n",
    "                                                                                                                  '8 Ciliated_IsthUTJ1'\n",
    "                                                                                                                 '9 Ciliated_IsthUTJ2'])]\n",
    "\n",
    "VP2h24hIAIU_batchcorrected_Fibroblast = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['3 Fibroblast_Pdgfra+',\n",
    "                                                                                                                  '4 Fibroblast_Pdgfra-'])]\n",
    "\n",
    "VP2h24hIAIU_batchcorrected_Muscle = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['12 Muscle'])]\n",
    "\n",
    "VP2h24hIAIU_batchcorrected_Immune = VP2h24hIAIU_batchcorrected[VP2h24hIAIU_batchcorrected.obs['leiden'].isin(['7 Immune1','13 Immune2'])]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Epithelial"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Secretory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Ciliated"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Fibroblast"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Muscle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "VP2h24hIAIU_batchcorrected_Immune"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Epithelial, svd_solver='arpack')\n",
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Secretory, svd_solver='arpack')\n",
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Ciliated, svd_solver='arpack')\n",
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Immune, svd_solver='arpack')\n",
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Fibroblast, svd_solver='arpack')\n",
    "sc.tl.pca(VP2h24hIAIU_batchcorrected_Muscle, svd_solver='arpack')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Epithelial, n_pcs=40, log=True)\n",
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Secretory, n_pcs=40, log=True)\n",
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Ciliated, n_pcs=40, log=True)\n",
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Immune, n_pcs=40, log=True)\n",
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Fibroblast, n_pcs=40, log=True)\n",
    "sc.pl.pca_variance_ratio(VP2h24hIAIU_batchcorrected_Muscle, n_pcs=40, log=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Epithelial, n_neighbors=10, n_pcs=40)\n",
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Secretory, n_neighbors=10, n_pcs=40)\n",
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Ciliated, n_neighbors=10, n_pcs=40)\n",
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Immune, n_neighbors=10, n_pcs=40)\n",
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Muscle, n_neighbors=10, n_pcs=40)\n",
    "sc.pp.neighbors(VP2h24hIAIU_batchcorrected_Fibroblast, n_neighbors=10, n_pcs=40)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Epithelial)\n",
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Secretory)\n",
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Ciliated)\n",
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Immune)\n",
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Muscle)\n",
    "sc.tl.umap(VP2h24hIAIU_batchcorrected_Fibroblast)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gene_dict = ['Epcam', 'Ovgp1', 'Foxj1', 'Myh11', 'Cd52', 'Ptprc', 'Pecam1', \n",
    "                  \n",
    "                  'Pdgfra', 'Twist2', 'Igf1', 'Dcn', 'Serpina1e', 'Myb', 'Pax8', 'Pax2', 'Cxcl17', 'Id2', 'Pdxk',\n",
    "             'Ccdc153', 'Il22ra1', 'Il22ra2', 'Stat3']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Epithelial, color=gene_dict)\n",
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Secretory, color=gene_dict)\n",
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Ciliated, color=gene_dict)\n",
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Fibroblast, color=gene_dict)\n",
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Muscle, color=gene_dict)\n",
    "sc.pl.umap(VP2h24hIAIU_batchcorrected_Immune, color=gene_dict)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Immune, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_batchcorrected_ImmuneCells_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_AllPops_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Epithelial, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_EpithelialOnly_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Secretory, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_SecretoryOnly_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Ciliated, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_CiliatedOnly_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Fibroblast, ['Il1b', 'Il10', 'Il15', 'Il18', 'Ifng', 'Tlr2', 'Tlr4'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_FibroblastOnly_immune genes.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected, ['Il10ra', 'Il10rb', 'Stat1', 'Stat3', 'Map2k1', 'Mapk8', 'Mapk14', 'Jak1', 'Tyk2'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_AllPops_IL22SignalingPathway.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pl.dotplot(VP2h24hIAIU_batchcorrected_Immune, ['Il10ra', 'Il10rb', 'Stat1', 'Stat3', 'Map2k1', 'Mapk8', 'Mapk14', 'Jak1', 'Tyk2'],\n",
    "groupby='Treatment', color_map='Blues', \n",
    "             save='_VP2h24hIAIU_batchcorrected_ImmuneCells_IL22SignalingPathway.pdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "winuthayanon-lab-191",
   "language": "python",
   "name": "winuthayanon-lab-191"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
